﻿using System;
using System.Collections.Generic;
using System.Drawing;
using System.IO;
using System.Runtime.Serialization;
using System.Runtime.Serialization.Formatters.Binary;
using MentalAlchemy;
using MentalAlchemy.Atomics;
using MentalAlchemy.Atomics.Types;

namespace MentalAlchemy.Molecules.MachineLearning
{
	/// <summary>
	/// Class for 2D Kohonen map.
	/// </summary>
	public class KohonenMap
	{
		#region - Variables. -
		public const int BATCH_TRAINING_EPOCHS = 100;

		protected List<List<Node2D>> nodes = new List<List<Node2D>>();
		protected Size size = new Size();
		protected List<int>[,] batches; 
		#endregion

		#region - Properties. -
		/// <summary>
		/// Defines size of the SOM.
		/// </summary>
		public Size Dimensions { get; set; }
		/// <summary>
		/// Training data.
		/// </summary>
		public List<float[]> TrainingData { get; set; }
		/// <summary>
		/// Array of training errors for each training epoch.
		/// </summary>
		public List<float> TrainingError { get; set; } 
		#endregion

		#region - Destructor. -
		~KohonenMap()
		{
			nodes.Clear();
			if (batches != null) Array.Clear(batches, 0, batches.Length);
			size = Size.Empty;

			if (TrainingData != null) TrainingData.Clear();
			if (TrainingError != null) TrainingError.Clear();
		}
		#endregion

		#region - Public methods. -
		/// <summary>
		/// Trains SOM using standard algorithm.
		/// </summary>
		/// <param name="epochs">Number of epochs to train.</param>
		/// <param name="earlyStop">Indicates whether early training stop is enabled.</param>
		public void Train(int epochs, bool earlyStop)
		{
			if (TrainingData != null)
			{
				// make initialization.
				InitNetwork();

				float sigma = 0.5f * Dimensions.Height, eta = 0.1f;	// see Haykin for parameters setting (page 585 in Russian translation).
				var t1 = (float)(Math.Log(sigma) / 1000.0);

				//WriteResultingWeights("epoch_0.txt");
				TrainingError = new List<float>();
				for (int i = 0; i < epochs; ++i)
				{
					float er = 0.0f;

					// run through the traning data.
					sigma *= (float)Math.Exp(-i * t1);
					eta *= (float)Math.Exp(-0.001 * i);
					float factor = CalculateFactor(sigma);
					foreach (var vector in TrainingData)
					{
						// get the node with the weights closest to the given vector.
						Node2D node = GetClosestNode(vector);

						var dists = StructMath.CalculateNeighborhoodDistancesExp(nodes, node, factor);
						CorrectWeights(vector, eta, dists);

						// calculate map error.
						er += CalculateError();
					}

					er /= TrainingData.Count;
					TrainingError.Add(er);

					// Debug output.
					//if (i > 0 && i % 10 == 0)
					//{
					//    WriteResultingWeights(String.Format("epoch_{0}.txt", i));
					//}

					// Stopping criterion: If error starts increasing then stop training.
					if (earlyStop && i > 10 && TrainingError[i] == TrainingError[i - 5])
					{
						break;
					}
				}

				// deactivate dead nodes.
				DeactivatePassiveNodes();
			}
		}

		/// <summary>
		/// Train SOM using batch algorithm.
		/// </summary>
		/// <param name="epochs">Number of training epochs.</param>
		/// <param name="earlyStop">Indicates whether early training stop is enabled.</param>
		public void TrainBatch(int epochs, bool earlyStop)
		{
			if (TrainingData != null)
			{
				// make initialization.
				InitNetwork();

				float sigma = 0.5f * Dimensions.Height, eta = 0.1f;	// see Haykin for parameters setting (page 585)
				var t1 = (float)(Math.Log(sigma) / 1000.0);

				//WriteResultingWeights("epoch_0.txt");
				TrainingError = new List<float>();
				int countValidBreak = 0;
				for (int i = 0; i < epochs; ++i)
				{
					float er = 0.0f;

					// find batches for the current nodes.
					batches = GetBatches();
					List<List<float[]>> avg = CalculateBatchCenters();

					// run through the traning data.
					sigma *= (float)Math.Exp(-i * t1);
					eta *= (float)Math.Exp(-0.001 * i);

					CorrectWeightsBatch(avg, sigma);

					er = CalculateError();
					if (float.IsNaN(er))
					{
						TrainingError.Add(100.0f);
					}
					TrainingError.Add(er);

					if (er == 0.0f) break;

					//if (i > 0 && i % 10 == 0)
					//{
					//    WriteResultingWeights(String.Format("epoch_{0}.txt", i));
					//}
					if (earlyStop && i > 10 && TrainingError[i] >= TrainingError[i - 1])
					{
						++countValidBreak;
						if (countValidBreak > 5)
						{
							break;
						}
					}

					// clear memory.
					avg.Clear();
				}

				// deactivate dead nodes.
				DeactivatePassiveNodes();
				GC.Collect();
			}
		}

		/// <summary>
		/// Create network of nodes spatially distributed on the 2d unit square and with random weights.
		/// </summary>
		public void InitNetwork()
		{
			if (Dimensions.IsEmpty) throw new Exception("[InitNetwork]: Dimensions.IsEmpty");
			if (TrainingData.Count == 0) throw new Exception("[InitNetwork]: TrainingData.Count == 0");

			int sizeTrainEntry = TrainingData[0].Length;
			float min = VectorMath.Min(TrainingData),
				  max = VectorMath.Max(TrainingData);

			float stepX = 1.0f / (Dimensions.Width - 1);
			float stepY = 1.0f / (Dimensions.Height - 1);

			var rnd = new Random();

			nodes.Clear();
			for (int i = 0; i < Dimensions.Height; ++i)
			{
				float row = i * stepY;

				nodes.Add(new List<Node2D>());
				nodes[i].Add(new Node2D());
				nodes[i][0].Coord = new PointF(0.0f, row);
				nodes[i][0].W = new List<float>(VectorMath.CreateRandomVector(rnd, sizeTrainEntry, min, max));
				nodes[i][0].IsActive = true;

				for (int j = 1; j < Dimensions.Width; ++j)
				{
					nodes[i].Add(new Node2D());
					nodes[i][j].Coord = new PointF((j * stepX), row);
					nodes[i][j].W = new List<float>(VectorMath.CreateRandomVector(rnd, sizeTrainEntry, min, max));
					nodes[i][j].IsActive = true;
				}
			}
		}

		/// <summary>
		/// Returns node, which is closest to the given vector using Euclidian measure.
		/// </summary>
		/// <param name="v">Vector.</param>
		/// <returns>Node closest to the vector [v].</returns>
		public Node2D GetClosestNode(float[] v)
		{
			var coords = GetClosestNodeCoords(v);
			return nodes[coords.Y][coords.X];
		}

		/// <summary>
		/// Gets coordinates of the node which is the closest to the given vector.
		/// </summary>
		/// <param name="v">Given vector.</param>
		/// <returns>Coordinates of the closest node.</returns>
		public Point GetClosestNodeCoords(float[] v)
		{
			var dist = StructMath.CalculateDistances(nodes, v);
			//var dist = (float[,])EvoPerformance.StructMath.CalculateDistances(nodes, v);
			return Performance.MatrixMath.IndexOfMin(dist);
		}

		/// <summary>
		/// Corrects weights using standard algorithm for SOM training.
		/// </summary>
		/// <param name="vector">Correcting vector.</param>
		/// <param name="eta">Learning speed.</param>
		/// <param name="dists">2D array of distances between SOM nodes and the given vector.</param>
		public void CorrectWeights(float[] vector, float eta, float[,] dists)
		{
			int height = dists.GetLength(0), width = dists.GetLength(1);
			for (int i = 0; i < height; i++)
			{
				for (int j = 0; j < width; j++)
				{
					float factor = eta * dists[i, j];
					Node2D node = nodes[i][j];
					for (int k = 0; k < vector.Length; ++k)
					{
						node.W[k] += factor * (vector[i] - node.W[i]);
					}
				}
			}
		}

		/// <summary>
		/// Corrects weights of all SOM nodes using information about batches.
		/// </summary>
		/// <param name="avgs">Averaged centers of batches.</param>
		/// <param name="sigma">Sigma parameter.</param>
		public void CorrectWeightsBatch(List<List<float[]>> avgs, float sigma)
		{
			int height = nodes.Count, width = nodes[0].Count;
			for (int i = 0; i < height; ++i)
			{
				for (int j = 0; j < width; ++j)
				{
					CorrectNodeBatch(nodes[i][j], avgs, sigma);
				}
			}
		}

		/// <summary>
		/// Calculates weighted sum distance to the previously calculated batches.
		/// </summary>
		/// <param name="dists">Distances to different batches.</param>
		/// <returns>Sum distance to all batches.</returns>
		public float CalculateSumDistanceBatch(float[,] dists)
		{
			int height = nodes.Count, width = nodes[0].Count;
			float dist = 0.0f;
			for (int i = 0; i < height; ++i)
			{
				for (int j = 0; j < width; ++j)
				{
					dist += dists[i, j] * batches[i, j].Count;
				}
			}
			return dist;
		} 

/// <summary>
        /// Corrects nodes weights using information about batches.
        /// </summary>
        /// <param name="node">Node to correct.</param>
        /// <param name="avgs">Array of averaged batches centers.</param>
        /// <param name="sigma">Neighborhood parameter.</param>
		public void CorrectNodeBatch(Node2D node, List<List<float[]>> avgs, float sigma)
		{
			float factor = CalculateFactor(sigma);
			var dists = StructMath.CalculateNeighborhoodDistancesExp(nodes, node, factor);

			//
			// calculate normalized distances to the other nodes.
			float sumDist = CalculateSumDistanceBatch(dists);
			if (sumDist == 0) {return;}

			float sumDist_1 = 1.0f/sumDist;
			if (float.IsInfinity(sumDist_1)) { return; }

			int height = nodes.Count, width = nodes[0].Count;
			var v = new float[node.W.Count];
			for (int i = 0; i < height; ++i)
			{
				for (int j = 0; j < width; ++j)
				{
					var tempV = (float[])avgs[i][j].Clone();
					float dist = batches[i,j].Count * dists[i,j] * sumDist_1;

					if (float.IsNaN(dist)) { return; }
					//if (float.IsNaN(dist)) { throw new Exception("CorrectNodeBatch: dist = NaN"); }

					VectorMath.Mul(ref tempV, dist);
					VectorMath.Accumulate(ref v, tempV);
				}
			}

			// update given node's weights.
			node.W = new List<float>(v);
		}

        /// <summary>
        /// Calculates error measure as mean square distance from
        ///     training sample to the closest node.
        /// </summary>
        /// <returns>Error value.</returns>
		public float CalculateError()
		{
        	return Performance.StructMath.CalculateError(nodes, TrainingData);

			// deprecated
			//float er = 0.0f;
			//foreach (var v in TrainingData)
			//{
			//    Node2D bum = GetClosestNode(v);

			//    er += VectorMath.EuclidianDistanceSqr(bum.W.ToArray(), v);
			//}

			//return er/TrainingData.Count;
		}

        /// <summary>
        /// Returns list of weights of the active nodes.
        /// </summary>
        /// <returns>List of active nodes weights</returns>
		public List<float[]> GetResultingWeights()
		{
			var res = new List<float[]>();
			foreach (List<Node2D> tNodes in nodes)
			{
				foreach (Node2D node in tNodes)
				{
					if (!node.IsActive) { continue; }

					res.Add(node.W.ToArray());
				}
			}	// foreach (List<Node2D> tNodes in nodes)

			return res;
		}

        /// <summary>
        /// Writes weights of active nodes into file with the given filename.
        /// </summary>
        /// <param name="filename">Name of the file.</param>
		public void WriteResultingWeights (string filename)
		{
			using (var writer = new StreamWriter(filename))
			{
				foreach (List<Node2D> tNodes in nodes)
				{
					foreach (Node2D node in tNodes)
					{
						if (!node.IsActive) {continue;}

						foreach (float w in node.W)
						{
							writer.Write("\t{0}", w);
						}
						//writer.Write("\t{0:D}", node.countCor);
						writer.Write("\t{0}\t{1}\n", node.Coord.X, node.Coord.Y);	// node coordinates
						//writer.Write("\t{0:D}\n", node.countCor);
					}
				}	// foreach (List<Node2D> tNodes in nodes)
			}	// using (StreamWriter writer = new StreamWriter(filename))
		}

        /// <summary>
        /// Calculates batches information -- each node is assigned an 1d array of training
        ///     data indices, for which this node is the closest one.
        /// </summary>
        /// <returns>2d array of lists containing training samples indices. Each cell in 2d array corresponds to the SOM node.</returns>
		public List<int>[,] GetBatches()
		{
			return Performance.StructMath.GetBatches (nodes, TrainingData);

			//if (nodes == null || nodes.Count == 0) throw new Exception("[GetBatches]: (nodes == null || nodes.Count == 0)");

			//#region - Create resulting data structure. -
			//var res = new List<List<List<int>>>();
			//res.AddRange(new List<List<int>>[nodes.Count]);
			//for (int i = 0; i < res.Count; ++i)
			//{
			//    res[i] = new List<List<int>>();
			//    res[i].AddRange(new List<int>[nodes[0].Count]);
			//    for (int j = 0; j < res[i].Count; ++j)
			//    {
			//        res[i][j] = new List<int>();
			//    }
			//} 
			//#endregion

			//for (int i = 0; i < TrainingData.Count; ++i )
			//{
			//    Point idx = GetClosestNodeCoords(TrainingData[i]);

			//    res[idx.Y][idx.X].Add(i);
			//}

			//return res;
		}

        /// <summary>
        /// Calculates batches centers for each SOM node. This method requires
        ///     GetBatches method to be called previously.
        /// </summary>
        /// <returns>2d array of batches centers.</returns>
		public List<List<float[]>> CalculateBatchCenters()
		{
			if (nodes == null || nodes.Count == 0) throw new Exception("[CalculateBatchCenters]: (nodes == null || nodes.Count == 0)");

            #region - Calculate resulting data structure. -
            var res = new List<List<float[]>>();
            int width = nodes[0].Count, height = nodes.Count;
            for (int i = 0; i < height; ++i)
            {
                res.Add(new List<float[]>());
                for (int j = 0; j < width; ++j)
                {
                    res[i].Add(new float[width]);
                }
            } 
            #endregion

            #region - Calculate centers of the batches for each SOM node. -
            int vectorSize = TrainingData[0].Length;
            for (int i = 0; i < height; ++i)
            {
                for (int j = 0; j < width; ++j)
                {
                    float[] temp = VectorMath.Zeros(vectorSize);
                    for (int k = 0; k < batches[i,j].Count; ++k)
                    {
                        float[] v = TrainingData[batches[i,j][k]];
                        VectorMath.Accumulate(ref temp, v);
                    }
                    VectorMath.Divide(ref temp, batches[i,j].Count);
                    res[i][j] = temp;
                }
            } 
            #endregion

			return res;
		}

        /// <summary>
        /// Sets all nodes, which are not closest to any training sample, deactivated.
        /// </summary>
		public void DeactivatePassiveNodes ()
		{
			var batch = GetBatches();

			int height = batch.GetLength(0), width = batch.GetLength(1);
			for (int i = 0; i < height; ++i )   // cycle through all nodes.
			{
				for (int j = 0; j < width; ++j)
				{
					if (batch[i,j].Count == 0)
					{
						nodes[i][j].IsActive = false;
					}
				}
			}
		}

		/// <summary>
		/// Saves Kohonen network into the specified file.
		/// </summary>
		/// <param name="filename">File name.</param>
		public void Save (string filename)
		{
			IFormatter formatter = new BinaryFormatter();
			Stream stream = new FileStream(filename, FileMode.Create, FileAccess.Write, FileShare.None);
			formatter.Serialize(stream, this);
			stream.Close();
		}
		#endregion

		#region - Static methods. -
		/// <summary>
		/// Calculates factor for normal distribution using given sigma value.
		///     See usages of this function to clarify its goal.
		/// </summary>
		/// <param name="sigma">Sigma value.</param>
		/// <returns>Factor value.</returns>
		private static float CalculateFactor(float sigma)
		{
			return 0.5f / (sigma * sigma);
		} 

		/// <summary>
		/// Loads Kohonen map from the specified file.
		/// </summary>
		/// <param name="filename">Filename.</param>
		/// <returns>Loaded Kohonen network.</returns>
		public static KohonenMap Load (string filename)
		{
			IFormatter formatter = new BinaryFormatter();
			Stream stream = new FileStream(filename, FileMode.Open, FileAccess.Read, FileShare.None);
			var kMap = (KohonenMap)formatter.Deserialize(stream);
			stream.Close();

			return kMap;
		}
		#endregion
	}
}